1. Описание датасета

В качестве сырых данных возьмем анонимизированные данные о поездках желтых такси в Нью-Йорке.

Сырые данные о поездках в жёлтом такси можно найти на сайте TLC: https://www1.nyc.gov/site/tlc/about/tlc-trip-record-data.page

Эти данные разбиты на файлы по месяцам.

Будем использовать информацию о поездках до мая 2016 года включительно, т.к они имеют более удобный вид для обработки.

Задача проекта — научиться предсказывать количество поездок на ближайшие сутки (24 часа) в каждом районе Нью-Йорка;

Для простоты мы поделим область Нью-Йорка сеткой 50 на 50. Поскольку в проекте нужно предсказывать 1 временной ряд, то мы выберем одну из центральных областей (под номером 1231), и будем предсказывать количество поездок именно для нее. Но в перспективе нужно работать со всеми ячейками.

Для того, чтобы это сделать, сырые данные необходимо агрегировать по часам и районам.

Агрегированные данные будут представлять собой почасовые временные ряды с количеством поездок из каждого района.

За основу проекта брался мой недоделанный проект для Courser'ы, поэтому возможно у кого то будет такая же идея.

2.Первичный анализ данных.

Как видно, файл очень большой, но нас интересуют только данные о месте начала поездки и о времени, поэтому можно считывать не все данные.

В каждом из файлов содержится следующая информация о поездках:

1. время начала поездки
2. время окончания поездки
5,6. долгота и широта точки начала поездки
9,10. долгота и широта точки окончания поездки
3. количество пассажиров
4. расстояние по счётчику
7. тип тарифа (одна из шести категорий)
11. способ оплаты (одна из шести категорий)
9. стоимость поездки по счётчику
12. налог на счётчик
13. доплата за поездки в пиковые часы и ночью
14. доплата за проезд по платным дорогам
17. доплата за проезд, взимаемая с каждой поездки с января 2015
15. размер чаевых
16. общая стоимость поездки
18. Общая сумма, взимаемая с пассажиров. Не включает в себя наличные чаевые.
9. бинарный флаг, показывающий, были ли данные о поездке получены немедленно после её окончания, или какое-то время хранились в памяти автомобиля.

За один месяц имеется почти 12 миллионов записей о поездках.

В основном данные типа float, т.к тут много информации о ценах, расстоянии и координатах.

Считаем тоже самое, но теперь только нужные столбцы.

Тем самым, мы сократили количество использованной памяти почти в 4 раза.

Предобработка данных

Теперь необходимо сделать агрегацию даннных по областям Нью-Йорка.

В данный момент данные имеют вид:

Удалим записи которые имеют некоректные данные

Сделаем нормальные индексы

Добавим нулевые столбцы, в которых в этом месяце не было поездок.

Теперь можно сохранить данные за месяц.

Чтобы было удобно обрабатывать файлы разных месяцев, оьбъединим все в 1 общую функцию.

Обработаем еще несколько месяцев

Визуальный анализ данных

Подключим нужные библиотеки

Будем использовать 4 месяца. С 2016-02 по 2016-05.

Для начала отрисуем график числа поездок из выбранной нами области.

Вывод : из этого графика явно видно выраженную дневную и недельную сезонность.

Это логично, люди чаще пользуются такси в вечерное время (как и в будние дни).

Поскольку ряд имеет сезонность, он не является стационарным.

Визуально в течении первых 3х месяцев среднее значение и дисперсия ряда выглядят постоянно.

Но вот в мае, значения становятся чуть меньшье (также как и дисперсия).

Возможно это связано с тем, что в мае уже достаточно тепло, чтобы ходить пешком).

Скорее всего это следствие годовой сезонности, которую мы не можем рассмотреть на отрезке в несколько месяцев.

Появляется отличная возможнось с помощью статистических критериев проверить гипотезу о равенстве среднего числа поездок за май месяц и числа поездок за предыдущие месяцы.

Посмотрим какое распределение имеет число поездок за весь этот период.

Видно, что оно имеет несколько пиков, а точнее 2 основных пика. Основная масса приходится на большие значения (около 1000).

Распределение совсем не похоже на нормальное.

Среднее сильно отличается от мод.

Можем сказать, что чаще всего люди используют такси либо "мало", либо "много".

Теперь взгляним на картину в общем и поиграемся с картами.

Отрисуем суммарное количество поездок за май.

Отрисуем среднее в час количество поездок за май.

Начиная с 5 поездок среднее будет показываться желтым цветом , c 50 - зеленым.

Получили красивую карту.

Судя по карте основная масса поездок приходится на Манхэттен. Что логично, т.к это деловой центр города.

Но также достаточное количество поездок приходится на аэропорты.

Ячеек из неадекватных мест(водоемов и т.д) вроде бы нет, так что скорее всего мы все сделали правильно.

Теперь разложим ряд на компоненты и отрисуем STL декомпозицию.

Вывод:

Визуально ряд не имеет ярко выраженного тренда. Значения колеблются в одних и тех же пределах.

Но в майский месяц значения более сжаты, будто дисперсия рядо уменьшилась.

Сезонная компонента состоит из дневной сезонности.

Недельняя сезонность ушла в тренд. Видно как тренд раскладывается на отрезки, где минимумы приходятся на выходные дни.

Проверка гипотез

  1. Проверим гипотеу о стационарности ряда с помощью Критерия Дики-Фулера.
  2. Убедимся в том что данные имеют не нормальное распределение.
  3. Проверим гипотезу о равенсте среднего числа поездок в мае и предыдущих месяцах.
  4. Попарно сравним средние разных месяцев (множественная проверка гипотез).

1. Стационарность ряда

Возьмем функцию с лекции.

Исходя из критерия Дики-Фулера, получается что ряд стационарен (т.к значение p-value 2.323501e-09).

Но мы знаем, что это не так, т.к есть сезонность.

Возможно, он не сработал из-за большой зашумленности данных.

2. Нормальность данных

Визуально, очевидно, что данные распределены не нормально, но давайте в этом убедимся с помощью критерия Шапиро-Уилка:

$H_0\colon$ данные распределны нормально.

$H_1\colon$ не нормально.

Критерий уверенно отвергает нулевую гипотезу в пользу альтернативу, с достигаемым уровнем значимости pvalue=1.4312167991741817e-31.

Мы убедились, что данные распределны не нормально.

3. Среднее число поездок

Сначала просто посмотрим на средние.

Видно что в мае среднее прилично меньше.

Дисперсии также отличаются.

Далее, будем проверять гипотезу о равенсте среднего числа поездок против двусторонней альтернативы. (Вместо того, чтобы брать одностороннюю альтернативу).

Это делается для того, чтобы мы не переобучалсь, т.к мы уже видели данные и знаем, что среднее скорее всего меньше.

Попробуем использовать Критерий Стьюдента, для проверки гипотезы о том, что среднее количество поездок за май не равно среднему количеству поездок за предыдущие месяцы.

H0: число заказов в мае равно числу заказов за 2-4 месяц.

H1: число заказов в мае не равно числу заказов за 2-4 месяц.

Хотя он и требует нормальности данных и равенства дисперсий, посмотрим на результаты.

Критерий уверенно отвергает нулевую гипотезу в пользу альтернативу, с достигаемым уровнем значимости pvalue=2.444402377784712e-05.

Критерий говорит, что средние значения все таки различны.

Но чтобы быть уверенными в результатах, используем критерий, который не требует никаких дполнительных условий от данных.

Воспользуем одним из непараметрических критериев - Критерием Манна-Уитни.

Критерий также уверенно отвергает нулевую гипотезу, с достигаемым уровнем значимости pvalue=1.2434458841687715e-07.

Можем сделать вывод: число поездок в мае значительно уменьшилось.

4. Множественная проверка гипотез

Проведем попарные тесты на равенство средних.

С помощью Критерия Манна-Уитни посчитали p-value, для гипотез о равенстве средних.

Теперь применим Метод Холма для поправки на множественную проверку гипотез.

Как видно, гипотезы о равенстве средних отклоняются, для случаев кроме 2-3, 2-4 месяца.

В остальных же случаях считаем, что средние не равны.

Интересно, что средние 2-3, 2-4 равны, а 3-4 нет.

Предсказание временного ряда.

Будем делать предсказания спроса на такси на ближайшие сутки.

Наивный baseline

Видим довольно большие ошибки. Попробуем что-нибудь посильнее.

Посчитаем качество тройного экспоненциального сглаживания на кросс-валидации и примем его за baseline.

Предсказания вообще выглядят довольно неадекватно, принимая значения даже меньше 0. Но посмотрим какое качество будет в среднем.

rmse меньше чем было раньше, посмотрим что можно с этим делать.

Линейная модель + Arima на остатки

Идея: давайте построим обычную линейную модель на придуманных нами регрессионных признаках. А информацию, которая останется в остатках от регресси, будем предсказывать с помощью ARIM'ы. Конечно, можно передать регрессионные признаки в Arima (c помощью параметра exog), но Arima и так учится слишком большое количество времени, а прирост качества должен из-за этого быть очень незначительный.

В рамках моделей ARIMA можно учесть только одну из сезонностей.

Обычно в таких случаях сезонность с самым маленьким периодом явно моделируют с помощью аримы, а все остальные учитывают за счёт регрессионной компоненты.

Недельную сезонность будем учитывать с помощью гармонических признаков.

Теперь сделаем небольшую кросс-валидацию, для того чтобы подобрать параметры

  1. K - для признаков $ sin([1,…,T] ∗ 2πi/168), \quad cos([1,…,T]∗2πi/168), i \in 1,...,K $, учитывающих недельную сезонность

  2. log_count - число лаговых фичей

Теперь график выглядит очень похоже. Параметры возьмем K = 12, lag_count = 8.

Обучим такую модель на 3 месяцах и посмотрим остатки.

Отрисуем график для трейна.

Анализ остатков

На самом деле остатки сильно напоминают белый шум. Явного тренда нет, сезонности тоже.

Проверим несмещенность остатков, с помощью критерия Стьюдента и стационарность с помощью критерия Дики-Фулера.

Как мы видим, остатки являются несмещенными и стационарными. Критерии уверенно это признают.

Вообще, на этом можно было бы остановиться, но давайте посмотрим, можно ли что-нибудь сделать с помощью Arim'ы

Отрисуем корролеллограммы.

Из графика автокорреляции все таки видно, что некоторая сезонность осталась в данных. Попробуем его продифференцировать.

Еще раз отрисуем королеллограмму.

Для сезонной производной:

Для обычной производной:

Остановимся на обычной проиводной, т.к королеллограмма выглядит более независимой, чем для сезонного случая.

Подбор параметров Arima

Из первого графика можно сделать вывод что:

Q = 1
q= 1 или 2

Из второго графика:

P = 1
p - тоже возьмем не большим (попробуем 0 - 2), хотя значимых значений намного больше.

И параметры дифференцирования:

d = 1
D = 0

Подберем параметры Arim'ы

В итоге лучшый набор параметров оказался

p = 1, q = 1, P = 1, Q = 0

Отрисуем финальные остатки.

Теперь это точно белый шум. Они также не смещены и стационарны.

Давайте совместим линейную модель и Arim'у и проверим качество на кросс-валидации.

Финальная модель.

Посмотрим на результаты.

Видно что качество не стало лучше.

Было - Score: 11.880454673188066 80.4179657335404

Посмотрим на скор глазами.

Видно, что в целом ошибка составляет примерно 10%, но например бывают дни, в которые модель ошибается на все 20%. (Не считая первые сплиты с маленьким трейном). Скорее всего, это является следствием того, что количество поездок на такси является очень шумным показателем, в том плане, что на него влияет очень много факторов, таких как, например, погода/пробки/праздники, которые не учитываются нашей моделью, но очень сильно могут влиять на спрос такси (значительно завышая/занижая его). Было бы интересно попробовать добавить подобные признаки, такие как предсказание синоптиков, в линейную модель).

Так что вполне возможно, что подобная точность в 11,5% является приемлимой.

Предсказание на один из последних дней выглядит как то так.

Вывод

В данной работе мы решали задачу предсказания числа поездок такси в одном из центральных районов Нью-Йорка.

Выводы которые мы сделали в ходе знакомства с данными.

1. Основня часть поездок на такси в Нью-Йорке осуществляется из Манхэтенна.

2. Число поездок обладает, как дневной, так и недельной сезонностью. Это хорошо видно из графика.
Также, мы можем предположить, что есть и годовая сезонность. 
Потому что число поездок в майе в среднем меньше (это мы проверили несколькими критериями), чем 
в предыдущие более холодные месяцы. Это кажется вполне логичным, т.к чем теплее, тем более вероятно человек выберет 
пройтись пешком, вместо того чтобы брать такси.

3. Поскольку ряд обладает сзонностью, то он не является стационарным. Однако, критерий Дики-Фулера, отвергает гипотезу     о нестационарности.

4. Также, мы показали, что распределение числа поездок двумодальное, т.е люди катаются либо мало, либо много.
И оно, очевидно, не нормальное.

5. С помощью множественной проверки гипотез, мы показали что средние для различных месяцев отличаются.
Средние равны, только для 2 и 3, 2 и 4 месяца.

Для предсказания ряда, мы использовали линейную модель.

Мы сделали некоторые регрессионные признаки (день недели, час дня, лаги, гармонические признаки).

Проверяли качество на кросс-валидации.

В итоге удалось получить качество MAPE = 11,88%

Остатки от линейной модели оказались стационарными, несмещенными и похожими на белый шум.

Также мы попробовали обучить SARIMAX на остатках от линейной модели. Это не дало значительного прироста в качестве (удалось достигнуть MAPE = 11,49%).

Если смотреть на качество на разных фолдах (в ходе кросс-валидации), то можно увидеть, что есть дни в которые качество значительно проседает. Скорее всего это следствие того, что целевая переменная - число поездок на такси, довольно зашумлена, т.к она зависит от многих случайных факторов (погоды, пробки и т.д). Поэтому, возможно, что для первого приближения это неплохое решение.

Варианты улучшения:

1. Добавить, например, погодные признаки, индикаторы праздников и т.д
2. Можно рассмотреть в качестве признаков число поездок в соседних ячейках, т.к скорее всего они сильно коррелируют.
Это могло бы дать нам еще больше информации.

В данной задаче, решение с помощью линейной модели мне кажется оптимальным в плане простоты построения модели и настройки параметров. Ведь, изначально задача состояла в предсказании ряда для каждой из ячеек, и вариант, например, с построением модели ARIMA(или любой другой модели в которой надо кропотливо настраивать параметры) для каждой из 2500 ячеек кажется очень сомнительным.

P.S. За возможные грамматические ошибки прошу прощения.